home *** CD-ROM | disk | FTP | other *** search
/ Visual Basic Source Code / Visual Basic Source Code.iso / vbsource / optivc16 / fitdemo.cpp < prev    next >
Encoding:
C/C++ Source or Header  |  1999-03-06  |  24.8 KB  |  481 lines
  1. /***********************   FITDEMO.CPP  **********************************
  2. *                                                                        *
  3. *            Data-fitting demo program for                               *
  4. *                  O p t i V e c                                         *
  5. *          with Borland C++ 4.x, 5.x, C++ Builder,                       *
  6. *           or Microsoft Visual C++ 5.0 or 6.0                           *
  7. *                                                                        *
  8. *       Copyright 1996-1999 by Martin Sander                             *
  9. *                                                                        *
  10. *                                                                        *
  11. *       This sample program is meant to provide you with some basic      *
  12. *       examples of code for the use of OptiVec's data-fitting routines. *
  13. *       Especially for the non-linear and multi-experiment functions,    *
  14. *       the user will have to adapt this code to his specific problem.   *
  15. *                                                                        *
  16. **************************************************************************/
  17. /*
  18. Borland C++, Command-line:
  19.     a) 32-bit: type
  20.                BCC32 -W fitdemo.cpp vcf3w.lib
  21.     b) 16-bit: type
  22.                BCC -ml -W fitdemo.cpp vcl3w.lib mcl3w.lib
  23. Borland C++, IDE:
  24.     1. Open the new-project menu with  Project/New.
  25.     2. Create a project FITDEMO in the OptiVec directory,
  26.        e.g., \bc\optivec.
  27.     3.a) 32-bit: Choose Application[EXE] for Win32 GUI, static linking,
  28.                  single-thread.
  29.     3.b) 16-bit: Choose Application[EXE] for Windows3.x, memory model LARGE.
  30.     4. Hit OK.
  31.     5. In the Project window that will now be on the screen, delete the
  32.        nodes FITDEMO.DEF and FITDEMO.RC.
  33.     6.a) 32-bit: Add the node VCF3W.LIB.
  34.     6.b) 16-bit: Add the nodes VCL3W.LIB and MCL3W.LIB.
  35.     7. In Options/Project, be sure the include-file search path and the
  36.        library search path both include the respective  OPTIVEC directories,
  37.        e.g.: \bc\optivec\include and bc\optivec\lib, resprectively
  38.     8. Compile and run.
  39. Microsoft Visual C++:
  40.     1. Create a new project as a "Win32 application".
  41.     2. In the project settings, C/C++, Code Generation,
  42.        verify that single-thread debug is chosen.
  43.     3. Add \OptiVec\include to the include-file search path.
  44.     4. Add the files  FITDEMO.CPP  and  OVVCSD.LIB to your project.
  45.     5. Compile and run.
  46. */
  47.  
  48. #include <windows.h>                    /* Compiler's include files */
  49. #include <stdio.h>
  50. #include <math.h>
  51. #include <string.h>
  52.  
  53. #include <VDstd.h>                      /* OptiVec include files */
  54. #include <VDmath.h>
  55. #include <MDstd.h>
  56. #include <VIstd.h>
  57. #include <Vgraph.h>
  58. HWND     hWndMain;
  59. int      vview;
  60. dVector  XExp, X2, YExp, YFit, YExp2, YFit2, YExp3, YFit3;
  61. ui       sizex;
  62. double   FitPars[7];   // the highest number of parameters we will have in our examples
  63. int      ParStatus[7];
  64. #define  polydeg 5     // refers to the polynomial fitting example
  65. #ifndef M_PI
  66.     #define M_PI  3.14159265358979323846
  67. #endif
  68. NEWMATHERR  // this macro has an effect only for 16-bit BC
  69.  
  70. LONG FAR PASCAL MainMessageHandler (HWND, UINT, WPARAM, LPARAM);
  71.  
  72.     // the following function is used with VD_linfit:
  73. static void PolyModel( dVector BasFuncs, double x, unsigned  nfuncs )
  74. {  /* This function fills the vector BasFuncs with powers of x.
  75.       VD_linfit will then determine the coefficient for each of these
  76.       powers.
  77.       Note that the coefficients do not occur in the model function!
  78.       The basis functions with known coefficients (whose fitting has
  79.       been disabled) are still part of the model and must be calculated.
  80.       You will see below that the values of the known (disabled)
  81.       coefficients  must be set prior to calling VD_linfit.          */
  82.      BasFuncs[0] = 1.0;
  83.      for( unsigned i=1; i<nfuncs; i++ )
  84.           BasFuncs[i] = BasFuncs[i-1]*x;
  85. }
  86.  
  87.    // the following function is used with VD_nonlinfit:
  88. static void VPolyModel( dVector Y, dVector X, ui size )
  89. {   /* Here, the model function has to fill a whole result vector,
  90.        using your first guess of FitPars. In contrast to the
  91.        linear case, now the coefficients are explicitly used in the
  92.        model function. You must initialize FitPars with something,
  93.        even if you have no idea about the result.
  94.        FitPars must be global, so that the model function can access the
  95.        parameters. With VD_nonlinfit, you can use just any functions
  96.        in your model.
  97.        For better comparison, we use the same polynomial approximation
  98.        as before, but now we code it as if we didn't know that a
  99.        polynomial actually is linear in its coefficients (and as if
  100.        we didn't know either how to efficiently code a polynomial at all).  */
  101.     double xi;
  102.     for( ui i=0; i<size; i++ )
  103.     {
  104.         xi  = X[i];
  105.         Y[i]= FitPars[0]
  106.             + FitPars[1] * xi
  107.             + FitPars[2] * xi * xi
  108.             + FitPars[3] * xi * xi * xi
  109.             + FitPars[4] * xi * xi * xi * xi
  110.             + FitPars[5] * xi * xi * xi * xi * xi;
  111.    }
  112. }
  113.  
  114.    // the following function is used with VD_multiNonlinfit:
  115. static void VSineModel( dVector Y, dVector X, ui size, unsigned nExperiment )
  116. {   /* According to the parameter nExperiment, the model function must
  117.        choose the correct parameters for the calculation of the model
  118.        function. The model function itself is the same for all experiments. */
  119.     double omega, phase, amp;
  120.     switch( nExperiment )
  121.     {
  122.         case 0: phase = FitPars[1];
  123.                 amp   = FitPars[4];
  124.                 break;
  125.         case 1: phase = FitPars[2];
  126.                 amp   = FitPars[5];
  127.                 break;
  128.         case 2: phase = FitPars[3];
  129.                 amp   = FitPars[6];
  130.     }
  131.     omega = FitPars[0];  // we assume this parameter to be the same
  132.     VDx_sin( Y, X, size, omega, phase, amp );
  133. }
  134.  
  135. int PASCAL WinMain(HINSTANCE hInstance, HINSTANCE hPrevInstance,
  136.                    LPSTR /* lpCmdLine */, int /* nCmdShow */)
  137. {
  138.     MSG       msg;                      /* MSG structure to pass to windows proc */
  139.     WNDCLASS  wc;
  140.     char      *AppName;
  141.  
  142.     AppName = "FitDemo";
  143.     if(!hPrevInstance)
  144.     {
  145.         wc.style      = CS_HREDRAW | CS_VREDRAW;
  146.         wc.lpfnWndProc= MainMessageHandler;
  147.         wc.cbClsExtra = 0;
  148.         wc.cbWndExtra = 0;
  149.         wc.hInstance  = hInstance;
  150.         wc.hIcon      = LoadIcon (hInstance, AppName);
  151.         wc.hCursor    = LoadCursor (NULL, IDC_ARROW);
  152.         wc.hbrBackground  = (HBRUSH) GetStockObject (WHITE_BRUSH);
  153.         wc.lpszMenuName   = AppName;
  154.         wc.lpszClassName  = AppName;
  155.         RegisterClass (&wc);
  156.     }
  157.                              /* create application's Main window:  */
  158.     hWndMain = CreateWindow (AppName,
  159.                              "OptiVec Data-Fitting Demo",
  160.                              WS_OVERLAPPEDWINDOW,
  161.                              CW_USEDEFAULT,     /* Use default X, Y, and width  */
  162.                              CW_USEDEFAULT,
  163.                              CW_USEDEFAULT,
  164.                              CW_USEDEFAULT,
  165.                              NULL,              /* Parent window's handle      */
  166.                              NULL,              /* Default to Class Menu       */
  167.                              hInstance,         /* Instance of window          */
  168.                              NULL);             /* Create struct for WM_CREATE */
  169.  
  170.  
  171.     if (hWndMain == NULL)
  172.     {
  173.         MessageBox(NULL, "Could not create window in WinMain", NULL, MB_ICONEXCLAMATION);
  174.         return (1);
  175.     }
  176.  
  177.     ShowWindow(hWndMain, SW_SHOWMAXIMIZED);     /* Display main window      */
  178.     UpdateWindow(hWndMain);
  179.  
  180.     while(GetMessage(&msg, NULL, 0, 0)) /* Main message loop */
  181.     {
  182.         TranslateMessage(&msg);
  183.         DispatchMessage(&msg);
  184.     }
  185.  
  186.     UnregisterClass (AppName, hInstance);
  187.     return (msg.wParam);
  188. }
  189.  
  190. LONG FAR PASCAL MainMessageHandler(HWND hWnd, UINT Message, WPARAM wParam, LPARAM lParam)
  191. {
  192.     HDC         vDC;
  193.     PAINTSTRUCT ps;
  194.     RECT        r;
  195.  
  196.     VD_NONLINFITOPTIONS Opt;
  197.     VD_EXPERIMENT       ExpList[3];
  198.     char                DataText[80];
  199.  
  200.  
  201.     switch (Message)                    /* Windows Message Loop           */
  202.     {
  203.         case WM_CREATE:
  204.             sizex   = 200;
  205.             vview   = 0;
  206.             XExp  = VD_vector( sizex );
  207.             YExp  = VD_vector( sizex );
  208.             YFit  = VD_vector( sizex );
  209.             YExp2 = VD_vector( sizex );
  210.             YExp3 = VD_vector( sizex );
  211.             YFit2 = VD_vector( sizex );
  212.             YFit3 = VD_vector( sizex );
  213.             break;
  214.  
  215.         case WM_RBUTTONDOWN:
  216.             if( ++vview > 5 ) vview = 0;
  217.             InvalidateRect( hWnd, NULL, TRUE ); // Demand new Paint to display the next view
  218.             break;
  219.  
  220.         case WM_PAINT:
  221.             vDC = BeginPaint(hWndMain, &ps);
  222.             V_initPlot( hWndMain, vDC );
  223.             GetClientRect( hWndMain, &r );
  224.             V_setPlotRegion( r.left+(r.right-r.left)/4, r.top, r.right, r.bottom );
  225.             switch( vview )
  226.             {
  227.                 default: break;
  228.                 case 0:
  229.                   TextOut( vDC, 10, 10,
  230.                   "This is a series of graphs illustrating OptiVec data-fitting functions.", 71 );
  231.                   TextOut( vDC, 10, 40,
  232.                   "Always press the right mouse button to see the next view!", 56 );
  233.                   TextOut( vDC, 10, 70,
  234.                   "Press [Alt] [F4] to end the demonstration.", 41 );
  235.                break;
  236.                case 1:
  237.                   TextOut( vDC, 2, 10, "Suppose these are", 17 );
  238.                   TextOut( vDC, 2, 30, "your experimental", 17 );
  239.                   TextOut( vDC, 2, 50, "data points.", 12 );
  240.                   TextOut( vDC, 2, 70, "(Actually, they consist", 23 );
  241.                   TextOut( vDC, 2, 90, "of a simple cubic with", 22 );
  242.                   TextOut( vDC, 2,110, "1% added noise.)", 16 );
  243.  
  244.                   VD_ramp( XExp, sizex, 0, 1.0/(sizex-1) ); // "experimental" x-axis from 0 to 1
  245.                   VD_cubic( YExp, XExp, sizex );        // fake "measured" y-data as y = x^3
  246.                   VD_noise( YFit, sizex, 1, 0.005 );
  247.                   VD_addV( YExp, YExp, YFit, sizex ); // add 1% peak-to-peak "experimental noise"
  248.                   VD_xyAutoPlot( XExp, YExp, sizex, PS_NULL+SY_CROSS, GREEN );
  249.                break;
  250.                case 2:   // fit your data to one of the simplest models, a polynomial
  251.                   TextOut( vDC, 2, 10, "The red curve is a", 18 );
  252.                   TextOut( vDC, 2, 30, "fifth-order polynomial", 22 );
  253.                   TextOut( vDC, 2, 50, "fitted to your data.", 20 );
  254.                   TextOut( vDC, 2, 70, "Without noise, the", 18 );
  255.                   TextOut( vDC, 2, 90, "coefficients should", 19 );
  256.                   TextOut( vDC, 2,110, "have been:", 10 );
  257.                   TextOut( vDC, 2,130, "{0, 0, 1.0, 0, 0}", 17 );
  258.  
  259.                   VD_polyfit( FitPars, polydeg, XExp, YExp, sizex );
  260.                   VD_poly( YFit, XExp, sizex, FitPars, polydeg ); // calculate fit curve
  261.                   TextOut( vDC, 2,160, "Actually, we got:", 17 );
  262.                   sprintf( DataText, "a0 = %+7.5lf", FitPars[0] );
  263.                   TextOut( vDC, 2,180, DataText, 13 );
  264.                   sprintf( DataText, "a1 = %+7.5lf", FitPars[1] );
  265.                   TextOut( vDC, 2,200, DataText, 13 );
  266.                   sprintf( DataText, "a2 = %+7.5lf", FitPars[2] );
  267.                   TextOut( vDC, 2,220, DataText, 13 );
  268.                   sprintf( DataText, "a3 = %+7.5lf", FitPars[3] );
  269.                   TextOut( vDC, 2,240, DataText, 13 );
  270.                   sprintf( DataText, "a4 = %+7.5lf", FitPars[4] );
  271.                   TextOut( vDC, 2,260, DataText, 13 );
  272.                   sprintf( DataText, "a5 = %+7.5lf", FitPars[5] );
  273.                   TextOut( vDC, 2,280, DataText, 13 );
  274.                   TextOut( vDC, 2,310, "Note how even moderate", 22 );
  275.                   TextOut( vDC, 2,330, "noise leads to rather", 21 );
  276.                   TextOut( vDC, 2,350, "large errors in the", 19 );
  277.                   TextOut( vDC, 2,370, "fit parameters, if", 18 );
  278.                   TextOut( vDC, 2,390, "there are too many", 18 );
  279.                   TextOut( vDC, 2,410, "'free' parameters.", 18 );
  280.                   VD_xy2AutoPlot( XExp, YExp, sizex, PS_NULL | SY_CROSS, GREEN,
  281.                                   XExp, YFit, sizex, PS_SOLID, LIGHTRED );
  282.                break;
  283.                case 3:  // refine your fit by switching to a general linear model,
  284.                         // giving you the chance to consider only the uneven terms
  285.                   TextOut( vDC, 2, 10, "Suppose you know that", 21 );
  286.                   TextOut( vDC, 2, 30, "the coefficients of", 19 );
  287.                   TextOut( vDC, 2, 50, "all even terms are 0.", 21 );
  288.                   TextOut( vDC, 2, 70, "Then you fit to your", 20 );
  289.                   TextOut( vDC, 2, 90, "own linear model,", 17 );
  290.                   TextOut( vDC, 2,110, "consisting only of", 18 );
  291.                   TextOut( vDC, 2,130, "uneven terms.", 13 );
  292.                   TextOut( vDC, 2,160, "Now we get:", 11 );
  293.  
  294.                   ParStatus[0] = ParStatus[2] = ParStatus[4] = 0;  //disable fitting of even terms
  295.                   FitPars[0] = FitPars[2] = FitPars[4] = 0.0;  // set their coefficients to the known value, 0.0
  296.                   ParStatus[1] = ParStatus[3] = ParStatus[5] = 1;  // enable fitting of uneven terms
  297.                      // the disabled fitting parameters must be initialized before calling VD_linfit !
  298.                   VD_linfit( FitPars, ParStatus, polydeg+1,
  299.                              XExp, YFit, sizex,
  300.                              PolyModel );
  301.                   VD_poly( YFit, XExp, sizex, FitPars, polydeg ); // calculate new fit curve
  302.  
  303.                   TextOut( vDC, 2,180, "a0 = 0 (fix)", 12 );
  304.                   sprintf( DataText, "a1 = %+7.5lf", FitPars[1] );
  305.                   TextOut( vDC, 2,200, DataText, 13 );
  306.                   TextOut( vDC, 2,220, "a2 = 0 (fix)", 12 );
  307.                   sprintf( DataText, "a3 = %+7.5lf", FitPars[3] );
  308.                   TextOut( vDC, 2,240, DataText, 13 );
  309.                   TextOut( vDC, 2,260, "a4 = 0 (fix)", 12 );
  310.                   sprintf( DataText, "a5 = %+7.5lf", FitPars[5] );
  311.                   TextOut( vDC, 2,280, DataText, 13 );
  312.                   TextOut( vDC, 2,310, "This is about as close", 22 );
  313.                   TextOut( vDC, 2,330, "as we can get in the", 20 );
  314.                   TextOut( vDC, 2,350, "presence of noise.", 18 );
  315.                   VD_xy2AutoPlot( XExp, YExp, sizex, PS_NULL | SY_CROSS, GREEN,
  316.                                   XExp, YFit, sizex, PS_SOLID, LIGHTRED );
  317.                   break;
  318.                case 4:  // here, we mis-use a non-linear fitting algorithm
  319.                         // for our simple problem.
  320.                   TextOut( vDC, 2, 10, "Let's fire with the", 19 );
  321.                   TextOut( vDC, 2, 30, "'cannon' of a non-", 18 );
  322.                   TextOut( vDC, 2, 50, "linear fit on our", 17 );
  323.                   TextOut( vDC, 2, 70, "simple 'sparrow'", 16 );
  324.                   TextOut( vDC, 2, 90, "problem.", 8 );
  325.                   TextOut( vDC, 2,110, "It takes much longer", 20 );
  326.                   TextOut( vDC, 2,130, "to find the result...", 21 );
  327.  
  328.                   ParStatus[0] = ParStatus[2] = ParStatus[4] = 0;  //disable fitting of even terms, as before
  329.                   ParStatus[1] = ParStatus[3] = ParStatus[5] = 1;  // enable fitting of uneven terms
  330.                   VD_getNonlinfitOptions( &Opt );
  331.                   Opt.FigureOfMerit = 0;  // choose least-square fitting
  332.                   Opt.AbsTolChi = 1.e-4;
  333.                   Opt.FracTolChi = 1.e-3;  // makes the fit fast, but not very accurate
  334.                   Opt.LevelOfMethod = 3;  // if you fear you might jump into a
  335.                      // local rather than the true global parameter optimum, try
  336.                      // LevelOfMethod = 7 - but only if you have time to wait for the result
  337.                   VD_setNonlinfitOptions( &Opt );
  338.  
  339.                   FitPars[0] = FitPars[2] = FitPars[4] = 0.0;  // set known coefficients to the value, 0.0, as before
  340.                   FitPars[1] = FitPars[3] = FitPars[5] = 1.5;  // you must provide some guess here!
  341.                         // all fitting parameters must be initialized before calling VD_nonlinfit !
  342.                   VD_nonlinfit( FitPars, ParStatus, polydeg+1,
  343.                                 XExp, YExp, sizex,
  344.                                 VPolyModel, NULL );
  345.                    // If you know the derivatives with respect to each parameter, put
  346.                    // your knowledge into a DerivModel function and replace the 'NULL'
  347.                    // parameter with it. (Actually, here we do know; but let's assume
  348.                    // we don't, and have VD_nonlinfit call the numeric differentiation
  349.                    // procedure.)
  350.                   VPolyModel( YFit, XExp, sizex );  // get fit curve from model
  351.  
  352.                   TextOut( vDC, 2,160, "But finally we get:", 19 );
  353.                   TextOut( vDC, 2,180, "a0 = 0 (fix)", 12 );
  354.                   sprintf( DataText, "a1 = %+7.5lf", FitPars[1] );
  355.                   TextOut( vDC, 2,200, DataText, 13 );
  356.                   TextOut( vDC, 2,220, "a2 = 0 (fix)", 12 );
  357.                   sprintf( DataText, "a3 = %+7.5lf", FitPars[3] );
  358.                   TextOut( vDC, 2,240, DataText, 13 );
  359.                   TextOut( vDC, 2,260, "a4 = 0 (fix)", 12 );
  360.                   sprintf( DataText, "a5 = %+7.5lf", FitPars[5] );
  361.                   TextOut( vDC, 2,280, DataText, 13 );
  362.                   TextOut( vDC, 2,310, "That is virtually the", 21 );
  363.                   TextOut( vDC, 2,330, "same as before.", 15 );
  364.                   VD_xy2AutoPlot( XExp, YExp, sizex, PS_NULL | SY_CROSS, GREEN,
  365.                                   XExp, YFit, sizex, PS_SOLID, LIGHTRED );
  366.                break;
  367.                case 5:  // finally, let's suppose you have several experimental
  368.                         // curves, measuring the same physical process under
  369.                         // slightly different conditions.
  370.                         // Say, we have a vibration, and each measurement
  371.                         // begins with a different phase and has a somewhat
  372.                         // different amplitude, but the same frequency.
  373.                   TextOut( vDC, 2, 10, "Finally, see how you", 20 );
  374.                   TextOut( vDC, 2, 30, "can fit several sets", 20 );
  375.                   TextOut( vDC, 2, 50, "of experimental data", 20 );
  376.                   TextOut( vDC, 2, 70, "at once.", 8 );
  377.                   TextOut( vDC, 2, 90, "This is much better", 19 );
  378.                   TextOut( vDC, 2,110, "than fitting them", 17 );
  379.                   TextOut( vDC, 2,130, "separately and averag-", 22 );
  380.                   TextOut( vDC, 2,150, "ing over the results.", 21 );
  381.                   TextOut( vDC, 2,170, "First you see the", 17 );
  382.                   TextOut( vDC, 2,190, "'experimental' data", 19 );
  383.                   TextOut( vDC, 2,210, "Please wait (may take", 21 );
  384.                   TextOut( vDC, 2,230, "several minutes)...", 19 );
  385.  
  386.                   VD_ramp( XExp, sizex, 0, 1.0/(sizex-1) ); // x-axis again from 0 to 1
  387.                   VDx_sin( YExp,  XExp, sizex, 15.0,  0.0, 1.2 );   // first "measurement"
  388.                   VDx_sin( YExp2, XExp, sizex, 15.0,  0.5, 1.0 );   // second "measurement"
  389.                   VDx_sin( YExp3, XExp, sizex, 15.0, -1.8, 0.75 ); // third "measurement"
  390.  
  391.                   VD_xy2AutoPlot( XExp, YExp,  sizex, PS_NULL + SY_CROSS, GREEN,
  392.                                   XExp, YExp2, sizex, PS_NULL + SY_CROSS, BLUE );
  393.                   VD_xyDataPlot(  XExp, YExp3, sizex, PS_NULL + SY_CROSS, MAGENTA );
  394.  
  395.                       // cram your experiments into the array of VD_EXPERIMENT structs:
  396.                   ExpList[0].X = XExp;   ExpList[0].Y = YExp;    ExpList[0].size = sizex;
  397.                   ExpList[1].X = XExp;   ExpList[1].Y = YExp2;   ExpList[1].size = sizex;
  398.                   ExpList[2].X = XExp;   ExpList[2].Y = YExp3;   ExpList[2].size = sizex;
  399.                        // we are not using the InvVar and WeightOfExperiment fields
  400.                        // of ExpList, as we are not weighting the data.
  401.  
  402.                   VI_equC( ParStatus, 7, 1 ); // we have 1 frequency, 3 phases, and 3 amplitudes,
  403.                                               // and all these 7 parameters are unknown.
  404.                                               // We must provide a first guess for each of them:
  405.                   FitPars[0] = 10.0;   // the frequency term
  406.                   FitPars[1] = FitPars[2] = FitPars[3] = 0.0;  // the three phases
  407.                   FitPars[4] = FitPars[5] = FitPars[6] = 1.5;  // the three amplitudes
  408.                   VD_getNonlinfitOptions( &Opt );
  409.                   Opt.AbsTolChi  = 1.e-8;
  410.                   Opt.FracTolChi = 1.e-6;  // force higher accuracy to avoid premature break-off
  411.                   /*  Unlike almost every other fitting routine available, you
  412.                       can get a result even for input parameters much farther off
  413.                       from the true value than the "guesses" chosen above.
  414.                       But then you must run VD_multiNonlinfit at "full power" and
  415.                       enable the following line:                               */
  416.                   // Opt.LevelOfMethod = 7;
  417.                   VD_setNonlinfitOptions( &Opt );
  418.  
  419.                   VD_multiNonlinfit( FitPars, ParStatus, 7,
  420.                                      ExpList, 3,
  421.                                      VSineModel,
  422.                                      NULL ); // Again, we pretend we don't know the derivatives
  423.                             // bring the phases into the range -PI < phase < + PI
  424.                   FitPars[1] = fmod( FitPars[1], 2.0*M_PI );
  425.                   if( FitPars[1] > M_PI ) FitPars[1] -= 2.0*M_PI;
  426.                   else if( FitPars[1] < -M_PI ) FitPars[1] += 2.0*M_PI;
  427.                   FitPars[2] = fmod( FitPars[2], 2.0*M_PI );
  428.                   if( FitPars[2] > M_PI ) FitPars[2] -= 2.0*M_PI;
  429.                   else if( FitPars[2] < -M_PI ) FitPars[2] += 2.0*M_PI;
  430.                   FitPars[3] = fmod( FitPars[3], 2.0*M_PI );
  431.                   if( FitPars[3] > M_PI ) FitPars[3] -= 2.0*M_PI;
  432.                   else if( FitPars[3] < -M_PI ) FitPars[3] += 2.0*M_PI;
  433.                   VSineModel( YFit,  XExp, sizex, 0 );  // get fit curves from your model
  434.                   VSineModel( YFit2, XExp, sizex, 1 );
  435.                   VSineModel( YFit3, XExp, sizex, 2 );
  436.                   #if defined __FLAT__ || defined _WIN32
  437.                      SetViewportOrgEx( vDC, r.left, r.top, NULL );
  438.                   #else
  439.                      SetViewportOrg( vDC, r.left, r.top );
  440.                   #endif  // to write text after plotting, go back to full window
  441.  
  442.                   TextOut( vDC, 2,265, "Here are the results", 20 );
  443.                   TextOut( vDC, 2,285, "(in brackets: 'true')", 21 );
  444.                   sprintf( DataText, "freq =%+7.5lf (15.0) ", FitPars[0] );
  445.                   TextOut( vDC, 2,320, DataText, 22 );
  446.                   sprintf( DataText, "ph1  = %+7.5lf (0.0)", FitPars[1] );
  447.                   TextOut( vDC, 2,340, DataText, 21 );
  448.                   sprintf( DataText, "ph2  = %+7.5lf (0.5)", FitPars[2] );
  449.                   TextOut( vDC, 2,360, DataText, 21 );
  450.                   sprintf( DataText, "ph3  = %+7.5lf (-1.8)", FitPars[3] );
  451.                   TextOut( vDC, 2,380, DataText, 22 );
  452.                   sprintf( DataText, "amp1 = %+7.5lf (1.2)", FitPars[4] );
  453.                   TextOut( vDC, 2,400, DataText, 21 );
  454.                   sprintf( DataText, "amp2 = %+7.5lf (1.0)", FitPars[5] );
  455.                   TextOut( vDC, 2,420, DataText, 21 );
  456.                   sprintf( DataText, "amp3 = %+7.5lf (0.75)", FitPars[6] );
  457.                   TextOut( vDC, 2,440, DataText, 22 );
  458.  
  459.                   V_continuePlot();   // go back to plot viewport
  460.                   VD_xyDataPlot(  XExp, YFit, sizex, PS_SOLID, LIGHTRED );
  461.                   VD_xyDataPlot(  XExp, YFit2, sizex, PS_SOLID, LIGHTRED );
  462.                   VD_xyDataPlot(  XExp, YFit3, sizex, PS_SOLID, LIGHTRED );
  463.             }
  464.             EndPaint(hWndMain, &ps);
  465.             break;
  466.  
  467.         case WM_CLOSE:
  468.             V_freeAll();                    /* delete all allocated vectors */
  469.             DestroyWindow(hWnd);
  470.             if (hWnd == hWndMain)
  471.                 PostQuitMessage(0);
  472.             break;
  473.  
  474.         default:
  475.             return (DefWindowProc(hWnd, Message, wParam, lParam));
  476.     }  // end of switch( message )
  477.     return (0);
  478. }
  479.  
  480.  
  481.